Magnet hypothesis: plant-pollinator interactions

Purpose: A test of the magnet hypothesis was examined in Mojave National Preserve by Ally Ruttan.

Hypothesis: Floral resource island created by shrubs and the associated beneficiary annual plants will positively and non-additively influence pollinator visitation rates.

Predictions:
(1) The frequency and duration of pollinator visitations to annuals is greater under shrubs than in the paired-open microsites (magnet H because of concentration).
(2) Annual plants under flowering entomophilous shrubs (Larrea tridentata) will have a higher frequency and duration of pollinator visitations than annual plants under anemophilous shrubs (Ambrosia dumosa) because of higher concentrations of suitable floral resources for pollinators (specificity of pollinator faciliation).
(3) Shrubs with annuals in their understory will have a higher frequency and duration of pollinator visitations than shrubs without annuals due to increased concentrations of floral resources for pollinators (reverse magnet effect and reciprocal benefits).
(4) Sites with both shrubs and annuals will have the highest frequency and duration of pollinator visitations to both the shrubs and the annuals (i.e. annuals under shrubs also with flowers are visited the most).

An interesting corollary is that there are appropriate floral resources for desert pollinators, that they discriminate, and that entomophilous and anemophilous shrubs facilitate flowering similarly.

Data wrangling

#libraries####
library(tidyverse)
library(DT)
library(lubridate)

#meta-data####
meta <- read_csv("data/meta-data.csv")
datatable(meta)

#error is SD

#data####
data.2015 <- read_csv("data/MNP.2015.csv")
data.2016 <- read_csv("data/MNP.2016.csv")

#merge
data <- rbind.data.frame(data.2015, data.2016)

#code treatment properly
data <- data %>% rename(net.treatment = treatment) #%>% na.omit(data) 

#keep key columns and minimize working dataframe
data <- data %>% select(-name, -plant, -start, -stop, -ID, -recorder) %>% filter(insect.RTU != "none")

#set year and rep as characters
data$year <- as.character(data$year)
data$rep <- as.character(data$rep)

#convert times to total seconds then to hour
data$total.duration <- (as.numeric(data$total.duration))/3600
data$visitation.duration <- (as.numeric(data$visitation.duration))/3600

#recode net.treatment column
data <- data %>% mutate(net.treatment = ifelse(net.treatment %in% c("SA"), "Larrea flowers with annuals", ifelse(net.treatment %in% c("SX"), "Larrea flowers without annuals", ifelse(net.treatment %in% c("SAA"), "Annual flowers under Larrea",ifelse(net.treatment %in% c("OA"), "Annual flowers in open",ifelse(net.treatment %in% c("AMB"), "Annual flowers under Ambrosia","NA"))))))

#frequency wrangled by RTU####
frequency <- data %>% group_by(year, day, net.treatment, rep, insect.RTU) %>% 
  summarise(net.time = sum(total.duration), mean.time = mean(total.duration), mean.temp = mean(temperature), mean.var.temp = mean(error), net.visitation = sum(visitation.duration), mean.visitation.duration = mean(visitation.duration), net.floral.density = sum(floral.density), mean.floral.density = mean(floral.density), insect.richness = n_distinct(insect.RTU), count = n())

#rates needed
frequency <- frequency %>% mutate(rate.per.flower = (count/mean.floral.density)) 

#view frequency data
datatable(frequency)

#separate net.treatment to individual vectors/factors
#Expt 1. Single-magnet hypothesis
expt1 <- frequency %>% filter(net.treatment == "Annual flowers under Larrea" | net.treatment == "Annual flowers in open" | net.treatment == "Annual flowers under Ambrosia")

expt1 <- expt1 %>% mutate(microsite = ifelse(net.treatment %in% c("Annual flowers under Larrea"), "Larrea", ifelse(net.treatment %in% c("Annual flowers under Ambrosia"), "Ambrosia", ifelse(net.treatment %in% c("Annual flowers in open"), "open", "NA"))))

#by years
expt1.2015 <- expt1 %>% filter(year == 2015)
expt1.2016 <- expt1 %>% filter(year == 2016)

#Expt2. Shrub flower-annual interactions through pollinators
expt2 <- frequency %>% filter(net.treatment == "Larrea flowers with annuals" | net.treatment == "Larrea flowers without annuals") 

expt2 <- expt2 %>% mutate(treatment = ifelse(net.treatment %in% c("Larrea flowers with annuals"), "Larrea + annuals", ifelse(net.treatment %in% c("Larrea flowers without annuals"), "Larrea - annuals", "NA")))

expt2.2015 <- expt2 %>% filter(year == 2015)
expt2.2016 <- expt2 %>% filter(year == 2016)

#Data structure: A total of 4 dataframes, two experiments, two years each

map <- read_csv("data/locations.csv")
map

## # A tibble: 5 × 3
##    site      lat      long
##   <int>    <dbl>     <dbl>
## 1     1 35.38890 -115.4050
## 2     2 35.31008 -115.3928
## 3     3 35.41113 -115.4016
## 4     4 35.41392 -115.3938
## 5     5 35.06128 -115.6644

Maps

#just one point
map <- map %>% filter(lat < 35.3)

require(ggmap)
cali <- get_map(location = c(lon = -115.6, lat = 35.06), zoom = 10)
#cali <-get_googlemap("california", crop= FALSE, zoom = 10)
p <-ggmap(cali)
p + geom_point(data=map, aes(x=long, y=lat), alpha = .25, size = 6, color = "blue")

#really zoomed out
cali <- get_map(location = c(lon = -115.6, lat = 35.06), zoom = 9)
p <-ggmap(cali)
p + geom_point(data=map, aes(x=long, y=lat), alpha = .25, size = 4, color = "blue")

Data visualization

#expt1####
#visitations####
#unweighted
ggplot(expt1, aes(microsite, rate.per.flower, fill = insect.RTU)) + geom_boxplot() + ylab("rate per flower") + scale_fill_brewer(palette = "Paired") + facet_grid(~year)

#weighted
ggplot(expt1, aes(microsite, rate.per.flower/mean.time, fill = insect.RTU)) + geom_boxplot() + ylab("rate per flower") + scale_fill_brewer(palette = "Paired") + facet_grid(~year)

#ggplot(expt1, aes(microsite, rate.per.flower/mean.time, fill = insect.RTU)) + geom_bar(stat="identity") + ylab("rate per flower") + scale_fill_brewer(palette = "Blues") + facet_grid(~year)

#visit durations####
#net visitation time
expt1 <- expt1 %>% filter(net.visitation < 2)
ggplot(expt1, aes(microsite, net.visitation, fill = insect.RTU)) + geom_boxplot() + ylab("net duration of visits per hour") + scale_fill_brewer(palette = "Paired") + facet_grid(~year)

expt1 <- expt1 %>% filter(net.visitation < 2)
ggplot(expt1, aes(microsite, net.visitation/mean.time, fill = insect.RTU)) + geom_boxplot() + ylab("net duration of visits per hour") + scale_fill_brewer(palette = "Paired") + facet_grid(~year)

#expt1 <- expt1 %>% filter(net.visitation < 2)
#ggplot(expt1, aes(microsite, net.visitation/mean.time, fill = insect.RTU)) + geom_bar(stat= "identity") + ylab("net duration of visits per hour") + scale_fill_brewer(palette = "Blues") + facet_grid(~year)

#temperature####
ggplot(expt1, aes(mean.temp, rate.per.flower/mean.time, color = insect.RTU)) + geom_point() + ylab("rate per flower") + geom_smooth(method = "lm") + facet_wrap(~year) + scale_color_brewer(palette = "Paired") 

#floral density####
ggplot(expt1, aes(net.floral.density, count, color = insect.RTU)) + geom_point() + ylab("count") + geom_smooth(method = "lm") + facet_wrap(~year) + scale_color_brewer(palette = "Paired")

#expt2####
#visitations####
expt2 <- expt2 %>% filter(rate.per.flower < 0.25)
ggplot(expt2, aes(treatment, rate.per.flower/mean.time, fill = insect.RTU)) + geom_boxplot() + ylab("rate per flower") + scale_fill_brewer(palette = "Paired") + facet_grid(~year)

#visit durations####
#net visitation time
expt2 <- expt2 %>% filter(net.visitation < 0.02)
ggplot(expt2, aes(treatment, net.visitation/mean.time, fill = insect.RTU)) + geom_boxplot() + ylab("net duration of visits per hour") + scale_fill_brewer(palette = "Paired") + facet_grid(~year)

#temperature####
ggplot(expt2, aes(mean.temp, rate.per.flower/mean.time, color = insect.RTU)) + geom_point() + ylab("rate per flower") + geom_smooth(method = "lm") + facet_wrap(~year) + scale_color_brewer(palette = "Paired")

#floral density####
ggplot(expt2, aes(net.floral.density, count, color = insect.RTU)) + geom_point() + ylab("count") + geom_smooth(method = "lm") + facet_wrap(~year) + scale_color_brewer(palette = "Paired")

Models

#Additive-term models for all insect taxa####
#expt1####
#2015####
#visitations
m <- glm(rate.per.flower~microsite*insect.RTU + mean.temp + offset(mean.time), family = "quasipoisson", data = expt1.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                      Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
## NULL                                   150     55.080              
## microsite             1   0.2363       149     54.843    0.3046    
## insect.RTU            2  10.9433       147     43.900 2.511e-11 ***
## mean.temp             1  15.2317       146     28.668 < 2.2e-16 ***
## microsite:insect.RTU  2   0.3105       144     28.358    0.5003    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

require(lsmeans)
lsmeans(m, pairwise~insect.RTU, adjust="tukey")

## $lsmeans
##  insect.RTU    lsmean        SE df asymp.LCL  asymp.UCL
##  bees       -1.068423 0.1025740 NA -1.269464 -0.8673812
##  flies      -1.555346 0.1349382 NA -1.819820 -1.2908717
##  other      -2.495249 0.2540097 NA -2.993099 -1.9973990
## 
## Results are averaged over the levels of: microsite 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast       estimate        SE df z.ratio p.value
##  bees - flies  0.4869232 0.1567628 NA   3.106  0.0054
##  bees - other  1.4268263 0.2682588 NA   5.319  <.0001
##  flies - other 0.9399031 0.2817638 NA   3.336  0.0025
## 
## Results are averaged over the levels of: microsite 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 3 estimates

#all contrasts
lsmeans(m, pairwise~microsite*insect.RTU, adjust="tukey")

## $lsmeans
##  microsite insect.RTU     lsmean        SE df  asymp.LCL  asymp.UCL
##  Larrea    bees       -0.2178995 0.1672474 NA -0.5456984  0.1098994
##  open      bees       -1.9189457 0.2128326 NA -2.3360901 -1.5018014
##  Larrea    flies      -0.5293760 0.2040048 NA -0.9292179 -0.1295340
##  open      flies      -2.5813157 0.2452641 NA -3.0620245 -2.1006069
##  Larrea    other      -1.4865626 0.3986665 NA -2.2679347 -0.7051906
##  open      other      -3.5039352 0.3639053 NA -4.2171764 -2.7906939
## 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                      estimate        SE df z.ratio p.value
##  Larrea,bees - open,bees      1.7010462 0.3231921 NA   5.263  <.0001
##  Larrea,bees - Larrea,flies   0.3114765 0.2347223 NA   1.327  0.7701
##  Larrea,bees - open,flies     2.3634162 0.3455851 NA   6.839  <.0001
##  Larrea,bees - Larrea,other   1.2686631 0.4116538 NA   3.082  0.0251
##  Larrea,bees - open,other     3.2860357 0.4362379 NA   7.533  <.0001
##  open,bees - Larrea,flies    -1.3895698 0.3402125 NA  -4.084  0.0006
##  open,bees - open,flies       0.6623700 0.2078334 NA   3.187  0.0180
##  open,bees - Larrea,other    -0.4323831 0.4888091 NA  -0.885  0.9503
##  open,bees - open,other       1.5849894 0.3438722 NA   4.609  0.0001
##  Larrea,flies - open,flies    2.0519397 0.3615405 NA   5.676  <.0001
##  Larrea,flies - Larrea,other  0.9571867 0.4294446 NA   2.229  0.2244
##  Larrea,flies - open,other    2.9745592 0.4491008 NA   6.623  <.0001
##  open,flies - Larrea,other   -1.0947531 0.5039085 NA  -2.173  0.2506
##  open,flies - open,other      0.9226195 0.3645236 NA   2.531  0.1151
##  Larrea,other - open,other    2.0173725 0.5697723 NA   3.541  0.0053
## 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 6 estimates

#require(multcomp)
#ph <- glht(m) #tests different from no effect
#summary(ph)

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -11.012, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##      alpha 
## -0.8063606

#rep at all as a nested, random effect
m <- glm(rate.per.flower~microsite*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "quasipoisson", data = expt1.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
## NULL                                       150     55.080              
## microsite                 1   0.2363       149     54.843  0.253977    
## mean.temp                 1  16.2764       148     38.567 < 2.2e-16 ***
## insect.RTU:rep           11  11.1940       137     27.373 4.566e-09 ***
## microsite:insect.RTU:rep 11   5.1218       126     22.251  0.003007 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lsmeans(m, pairwise~microsite*insect.RTU, adjust="tukey")

## $lsmeans
##  microsite insect.RTU     lsmean        SE df asymp.LCL   asymp.UCL
##  Larrea    bees       -0.3322604 0.1623242 NA -0.650410 -0.01411072
##  open      bees       -1.9352388 0.1958070 NA -2.319013 -1.55146412
##  Larrea    flies      -0.6743964 0.2031691 NA -1.072601 -0.27619220
##  open      flies      -2.5806671 0.2276444 NA -3.026842 -2.13449222
##  Larrea    other      -1.4521590 0.4045215 NA -2.245007 -0.65931137
##  open      other      -3.5216435 0.3544505 NA -4.216354 -2.82693325
## 
## Results are averaged over the levels of: rep 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                      estimate        SE df z.ratio p.value
##  Larrea,bees - open,bees      1.6029784 0.3013352 NA   5.320  <.0001
##  Larrea,bees - Larrea,flies   0.3421360 0.2346161 NA   1.458  0.6910
##  Larrea,bees - open,flies     2.2484067 0.3232723 NA   6.955  <.0001
##  Larrea,bees - Larrea,other   1.1198986 0.4148244 NA   2.700  0.0752
##  Larrea,bees - open,other     3.1893831 0.4204170 NA   7.586  <.0001
##  open,bees - Larrea,flies    -1.2608424 0.3220041 NA  -3.916  0.0013
##  open,bees - open,flies       0.6454282 0.1994691 NA   3.236  0.0154
##  open,bees - Larrea,other    -0.4830798 0.4860388 NA  -0.994  0.9201
##  open,bees - open,other       1.5864047 0.3414675 NA   4.646  <.0001
##  Larrea,flies - open,flies    1.9062707 0.3425952 NA   5.564  <.0001
##  Larrea,flies - Larrea,other  0.7777626 0.4340622 NA   1.792  0.4710
##  Larrea,flies - open,other    2.8472471 0.4355912 NA   6.537  <.0001
##  open,flies - Larrea,other   -1.1285081 0.5000046 NA  -2.257  0.2120
##  open,flies - open,other      0.9409764 0.3601134 NA   2.613  0.0941
##  Larrea,other - open,other    2.0694845 0.5672400 NA   3.648  0.0036
## 
## Results are averaged over the levels of: rep 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 6 estimates

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -14.459, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##      alpha 
## -0.8716176

#visit duration####
m <- glm(net.visitation~microsite*insect.RTU + mean.temp + offset(mean.time), family = "gaussian", data = expt1.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                      Df Deviance Resid. Df Resid. Dev Pr(>Chi)  
## NULL                                   150     34.958           
## microsite             1  0.44985       149     34.508  0.15488  
## insect.RTU            2  1.18290       147     33.325  0.06991 .
## mean.temp             1  0.85343       146     32.472  0.05007 .
## microsite:insect.RTU  2  0.45971       144     32.012  0.35560  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#lsmeans(m, pairwise~microsite*insect.RTU, adjust="tukey")

#rep at all as a nested, random effect
m <- glm(net.visitation~microsite*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "gaussian", data = expt1.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df Deviance Resid. Df Resid. Dev Pr(>Chi)  
## NULL                                       150     34.958           
## microsite                 1   0.4499       149     34.508  0.14618  
## mean.temp                 1   0.7730       148     33.735  0.05679 .
## insect.RTU:rep           11   2.4830       137     31.252  0.39010  
## microsite:insect.RTU:rep 11   4.4100       126     26.842  0.03661 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lsmeans(m, pairwise~microsite*insect.RTU, adjust="tukey")

## $lsmeans
##  microsite insect.RTU      lsmean        SE df   asymp.LCL asymp.UCL
##  Larrea    bees        0.25850843 0.1059353 NA  0.05087906 0.4661378
##  open      bees        0.28713323 0.1024138 NA  0.08640580 0.4878607
##  Larrea    flies      -0.03164143 0.1171209 NA -0.26119412 0.1979113
##  open      flies       0.28027390 0.1012677 NA  0.08179292 0.4787549
##  Larrea    other      -0.07348001 0.1463836 NA -0.36038660 0.2134266
##  open      other       0.12845277 0.1180315 NA -0.10288472 0.3597903
## 
## Results are averaged over the levels of: rep 
## Results are given on the identity (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                        estimate        SE df z.ratio p.value
##  Larrea,bees - open,bees     -0.028624793 0.1718683 NA  -0.167  1.0000
##  Larrea,bees - Larrea,flies   0.290149860 0.1342408 NA   2.161  0.2560
##  Larrea,bees - open,flies    -0.021765468 0.1703763 NA  -0.128  1.0000
##  Larrea,bees - Larrea,other   0.331988448 0.1460876 NA   2.273  0.2053
##  Larrea,bees - open,other     0.130055668 0.1816694 NA   0.716  0.9801
##  open,bees - Larrea,flies     0.318774652 0.1736033 NA   1.836  0.4423
##  open,bees - open,flies       0.006859325 0.1194568 NA   0.057  1.0000
##  open,bees - Larrea,other     0.360613241 0.2039921 NA   1.768  0.4868
##  open,bees - open,other       0.158680461 0.1330020 NA   1.193  0.8404
##  Larrea,flies - open,flies   -0.311915327 0.1723212 NA  -1.810  0.4591
##  Larrea,flies - Larrea,other  0.041838589 0.1630189 NA   0.257  0.9998
##  Larrea,flies - open,other   -0.160094192 0.1832970 NA  -0.873  0.9529
##  open,flies - Larrea,other    0.353753916 0.2025734 NA   1.746  0.5009
##  open,flies - open,other      0.151821136 0.1330203 NA   1.141  0.8641
##  Larrea,other - open,other   -0.201932780 0.2123279 NA  -0.951  0.9330
## 
## Results are averaged over the levels of: rep 
## P value adjustment: tukey method for comparing a family of 6 estimates

#2016####
#visitations
m <- glm(rate.per.flower~microsite*insect.RTU + mean.temp + offset(mean.time), family = "quasipoisson", data = expt1.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                      Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
## NULL                                   159    147.770              
## microsite             2    1.373       157    146.397    0.2137    
## insect.RTU            2   73.156       155     73.241 < 2.2e-16 ***
## mean.temp             1    7.848       154     65.393 2.662e-05 ***
## microsite:insect.RTU  4    1.007       150     64.386    0.6874    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lsmeans(m, pairwise~insect.RTU, adjust="tukey")

## $lsmeans
##  insect.RTU     lsmean         SE df  asymp.LCL  asymp.UCL
##  bees        0.2496914 0.06475311 NA  0.1227776  0.3766051
##  flies      -1.4164450 0.18215414 NA -1.7734606 -1.0594295
##  other      -2.0625151 0.42482659 NA -2.8951600 -1.2298703
## 
## Results are averaged over the levels of: microsite 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast       estimate        SE df z.ratio p.value
##  bees - flies  1.6661364 0.1923802 NA   8.661  <.0001
##  bees - other  2.3122065 0.4296771 NA   5.381  <.0001
##  flies - other 0.6460701 0.4621849 NA   1.398  0.3420
## 
## Results are averaged over the levels of: microsite 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 3 estimates

#all contrasts
lsmeans(m, pairwise~microsite*insect.RTU, adjust="tukey")

## $lsmeans
##  microsite insect.RTU     lsmean        SE df  asymp.LCL  asymp.UCL
##  Ambrosia  bees        0.4148953 0.1082462 NA  0.2027365  0.6270540
##  Larrea    bees        0.4734162 0.1048574 NA  0.2678995  0.6789328
##  open      bees       -0.1392373 0.1345258 NA -0.4029031  0.1244284
##  Ambrosia  flies      -1.2834767 0.3364496 NA -1.9429058 -0.6240476
##  Larrea    flies      -1.2602674 0.3046985 NA -1.8574655 -0.6630692
##  open      flies      -1.7055910 0.3110168 NA -2.3151726 -1.0960093
##  Ambrosia  other      -1.7337139 0.5848053 NA -2.8799112 -0.5875166
##  Larrea    other      -2.7066033 0.9365425 NA -4.5421928 -0.8710137
##  open      other      -1.7472283 0.6378914 NA -2.9974725 -0.4969841
## 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                           estimate        SE df z.ratio p.value
##  Ambrosia,bees - Larrea,bees     -0.05852090 0.1497825 NA  -0.391  1.0000
##  Ambrosia,bees - open,bees        0.55413258 0.1759813 NA   3.149  0.0434
##  Ambrosia,bees - Ambrosia,flies   1.69837194 0.3530132 NA   4.811  0.0001
##  Ambrosia,bees - Larrea,flies     1.67516261 0.3229430 NA   5.187  <.0001
##  Ambrosia,bees - open,flies       2.12048622 0.3311559 NA   6.403  <.0001
##  Ambrosia,bees - Ambrosia,other   2.14860913 0.5942802 NA   3.615  0.0091
##  Ambrosia,bees - Larrea,other     3.12149851 0.9429053 NA   3.311  0.0261
##  Ambrosia,bees - open,other       2.16212355 0.6473125 NA   3.340  0.0237
##  Larrea,bees - open,bees          0.61265348 0.1768830 NA   3.464  0.0156
##  Larrea,bees - Ambrosia,flies     1.75689284 0.3516082 NA   4.997  <.0001
##  Larrea,bees - Larrea,flies       1.73368351 0.3214501 NA   5.393  <.0001
##  Larrea,bees - open,flies         2.17900712 0.3317182 NA   6.569  <.0001
##  Larrea,bees - Ambrosia,other     2.20713004 0.5932583 NA   3.720  0.0062
##  Larrea,bees - Larrea,other       3.18001941 0.9426376 NA   3.374  0.0212
##  Larrea,bees - open,other         2.22064445 0.6470264 NA   3.432  0.0174
##  open,bees - Ambrosia,flies       1.14423937 0.3655781 NA   3.130  0.0459
##  open,bees - Larrea,flies         1.12103003 0.3362214 NA   3.334  0.0242
##  open,bees - open,flies           1.56635364 0.3243650 NA   4.829  <.0001
##  open,bees - Ambrosia,other       1.59447656 0.6036639 NA   2.641  0.1696
##  open,bees - Larrea,other         2.56736593 0.9451455 NA   2.716  0.1419
##  open,bees - open,other           1.60799097 0.6495473 NA   2.476  0.2438
##  Ambrosia,flies - Larrea,flies   -0.02320933 0.4533181 NA  -0.051  1.0000
##  Ambrosia,flies - open,flies      0.42211428 0.4608740 NA   0.916  0.9922
##  Ambrosia,flies - Ambrosia,other  0.45023719 0.6738581 NA   0.668  0.9991
##  Ambrosia,flies - Larrea,other    1.42312656 0.9953901 NA   1.430  0.8866
##  Ambrosia,flies - open,other      0.46375161 0.7217336 NA   0.643  0.9994
##  Larrea,flies - open,flies        0.44532361 0.4379378 NA   1.017  0.9844
##  Larrea,flies - Ambrosia,other    0.47344652 0.6586680 NA   0.719  0.9985
##  Larrea,flies - Larrea,other      1.44633589 0.9850853 NA   1.468  0.8704
##  Larrea,flies - open,other        0.48696094 0.7074314 NA   0.688  0.9989
##  open,flies - Ambrosia,other      0.02812292 0.6657857 NA   0.042  1.0000
##  open,flies - Larrea,other        1.00101229 0.9858165 NA   1.015  0.9845
##  open,flies - open,other          0.04163733 0.7073782 NA   0.059  1.0000
##  Ambrosia,other - Larrea,other    0.97288937 1.1045410 NA   0.881  0.9940
##  Ambrosia,other - open,other      0.01351442 0.8662347 NA   0.016  1.0000
##  Larrea,other - open,other       -0.95937495 1.1328599 NA  -0.847  0.9954
## 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 9 estimates

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -8.2587, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##     alpha 
## -0.627083

#require(multcomp)
#ph <- glht(m) #tests different from no effect
#summary(ph)

#rep at all as a nested, random effect
m <- glm(rate.per.flower~microsite*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "quasipoisson", data = expt1.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df Deviance Resid. Df Resid. Dev  Pr(>Chi)    
## NULL                                       158    147.684              
## microsite                 2    1.426       156    146.258    0.1991    
## mean.temp                 1    9.910       155    136.348 2.178e-06 ***
## insect.RTU:rep           12   76.826       143     59.522 < 2.2e-16 ***
## microsite:insect.RTU:rep 19    3.779       124     55.743    0.9802    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -10.524, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##      alpha 
## -0.7075865

#visit duration####
m <- glm(net.visitation~microsite*insect.RTU + mean.temp + offset(mean.time), family = "gaussian", data = expt1.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                      Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL                                   158     884.41         
## microsite             2   19.362       156     865.04   0.1694
## insect.RTU            2   23.123       154     841.92   0.1200
## mean.temp             1   12.559       153     829.36   0.1291
## microsite:insect.RTU  4   16.914       149     812.45   0.5409

#lsmeans(m, pairwise~microsite*insect.RTU, adjust="tukey")

#rep at all as a nested, random effect
m <- glm(net.visitation~microsite*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "gaussian", data = expt1.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df Deviance Resid. Df Resid. Dev Pr(>Chi)
## NULL                                       157     884.40         
## microsite                 2   19.458       155     864.94   0.1996
## mean.temp                 1   14.256       154     850.69   0.1244
## insect.RTU:rep           12   44.319       142     806.37   0.8344
## microsite:insect.RTU:rep 19   63.698       123     742.67   0.9381

#expt2####
#2015####
#visitations
m <- glm(rate.per.flower~net.treatment*insect.RTU + mean.temp + offset(mean.time), family = "quasipoisson", data = expt2.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df Deviance Resid. Df Resid. Dev Pr(>Chi)   
## NULL                                       124     2.5996            
## net.treatment             1 0.042087       123     2.5575 0.186673   
## insect.RTU                2 0.223468       121     2.3341 0.009763 **
## mean.temp                 1 0.000177       120     2.3339 0.931692   
## net.treatment:insect.RTU  2 0.049217       118     2.2847 0.360764   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lsmeans(m, pairwise~insect.RTU)

## $lsmeans
##  insect.RTU    lsmean        SE df asymp.LCL asymp.UCL
##  bees       -3.406125 0.1035943 NA -3.609166 -3.203084
##  flies      -3.858280 0.1604478 NA -4.172752 -3.543808
##  other      -4.115501 0.2666536 NA -4.638133 -3.592870
## 
## Results are averaged over the levels of: net.treatment 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast       estimate        SE df z.ratio p.value
##  bees - flies  0.4521551 0.1906980 NA   2.371  0.0466
##  bees - other  0.7093762 0.2852204 NA   2.487  0.0344
##  flies - other 0.2572211 0.3108413 NA   0.827  0.6859
## 
## Results are averaged over the levels of: net.treatment 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 3 estimates

#all contrasts
lsmeans(m, pairwise~net.treatment*insect.RTU, adjust="tukey")

## $lsmeans
##  net.treatment                  insect.RTU    lsmean        SE df
##  Larrea flowers with annuals    bees       -3.443655 0.1490332 NA
##  Larrea flowers without annuals bees       -3.368595 0.1431138 NA
##  Larrea flowers with annuals    flies      -4.065369 0.2634157 NA
##  Larrea flowers without annuals flies      -3.651191 0.1831933 NA
##  Larrea flowers with annuals    other      -4.505571 0.4485553 NA
##  Larrea flowers without annuals other      -3.725431 0.2878097 NA
##  asymp.LCL asymp.UCL
##  -3.735754 -3.151555
##  -3.649093 -3.088097
##  -4.581655 -3.549084
##  -4.010243 -3.292139
##  -5.384723 -3.626419
##  -4.289528 -3.161335
## 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                                                                   
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,bees     
##  Larrea flowers with annuals,bees - Larrea flowers with annuals,flies       
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,flies    
##  Larrea flowers with annuals,bees - Larrea flowers with annuals,other       
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,other    
##  Larrea flowers without annuals,bees - Larrea flowers with annuals,flies    
##  Larrea flowers without annuals,bees - Larrea flowers without annuals,flies 
##  Larrea flowers without annuals,bees - Larrea flowers with annuals,other    
##  Larrea flowers without annuals,bees - Larrea flowers without annuals,other 
##  Larrea flowers with annuals,flies - Larrea flowers without annuals,flies   
##  Larrea flowers with annuals,flies - Larrea flowers with annuals,other      
##  Larrea flowers with annuals,flies - Larrea flowers without annuals,other   
##  Larrea flowers without annuals,flies - Larrea flowers with annuals,other   
##  Larrea flowers without annuals,flies - Larrea flowers without annuals,other
##  Larrea flowers with annuals,other - Larrea flowers without annuals,other   
##     estimate        SE df z.ratio p.value
##  -0.07505962 0.2060529 NA  -0.364  0.9992
##   0.62171466 0.3023240 NA   2.056  0.3105
##   0.20753595 0.2360821 NA   0.879  0.9516
##   1.06191627 0.4725483 NA   2.247  0.2162
##   0.28177652 0.3226703 NA   0.873  0.9529
##   0.69677428 0.2994948 NA   2.326  0.1833
##   0.28259557 0.2324011 NA   1.216  0.8293
##   1.13697589 0.4707307 NA   2.415  0.1508
##   0.35683614 0.3201738 NA   1.115  0.8756
##  -0.41417871 0.3208135 NA  -1.291  0.7901
##   0.44020161 0.5201041 NA   0.846  0.9587
##  -0.33993814 0.3892815 NA  -0.873  0.9529
##   0.85438032 0.4845069 NA   1.763  0.4896
##   0.07424057 0.3409856 NA   0.218  0.9999
##  -0.78013975 0.5325937 NA  -1.465  0.6868
## 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 6 estimates

#require(multcomp)
#ph <- glht(m) #tests different from no effect
#summary(ph)

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -14.334, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##      alpha 
## -0.8182304

#rep at all as a nested, random effect
m <- glm(rate.per.flower~net.treatment*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "quasipoisson", data = expt2.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                              Df Deviance Resid. Df Resid. Dev  Pr(>Chi)
## NULL                                           124     2.5996          
## net.treatment                 1  0.04209       123     2.5575 0.1357642
## mean.temp                     1  0.00002       122     2.5575 0.9754898
## insect.RTU:rep               11  0.62872       111     1.9288 0.0004805
## net.treatment:insect.RTU:rep 11  0.15219       100     1.7766 0.7091117
##                                 
## NULL                            
## net.treatment                   
## mean.temp                       
## insect.RTU:rep               ***
## net.treatment:insect.RTU:rep    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lsmeans(m, pairwise~insect.RTU)

## $lsmeans
##  insect.RTU    lsmean         SE df asymp.LCL asymp.UCL
##  bees       -3.460593 0.09885027 NA -3.654336 -3.266851
##  flies      -3.937371 0.16788781 NA -4.266425 -3.608316
##  other      -4.447086 0.41054131 NA -5.251732 -3.642440
## 
## Results are averaged over the levels of: net.treatment, rep 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast       estimate        SE df z.ratio p.value
##  bees - flies  0.4767770 0.1948062 NA   2.447  0.0382
##  bees - other  0.9864926 0.4219121 NA   2.338  0.0507
##  flies - other 0.5097156 0.4435196 NA   1.149  0.4838
## 
## Results are averaged over the levels of: net.treatment, rep 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 3 estimates

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -17.432, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##      alpha 
## -0.8684939

#visit duration####
m <- glm(net.visitation~net.treatment*insect.RTU + mean.temp + offset(mean.time), family = "gaussian", data = expt2.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df Deviance Resid. Df Resid. Dev Pr(>Chi)   
## NULL                                       125     23.874            
## net.treatment             1  0.19369       124     23.680 0.304446   
## insect.RTU                2  0.23689       122     23.443 0.524709   
## mean.temp                 1  1.43728       121     22.006 0.005151 **
## net.treatment:insect.RTU  2  0.15015       119     21.856 0.664464   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#lsmeans(m, pairwise~net.treatment*insect.RTU, adjust="tukey")

#rep at all as a nested, random effect
m <- glm(net.visitation~net.treatment*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "gaussian", data = expt2.2015)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                              Df Deviance Resid. Df Resid. Dev Pr(>Chi)   
## NULL                                           125     23.874            
## net.treatment                 1  0.19369       124     23.680 0.315529   
## mean.temp                     1  1.40896       123     22.271 0.006789 **
## insect.RTU:rep               11  1.98162       112     20.289 0.503062   
## net.treatment:insect.RTU:rep 11  0.86994       101     19.419 0.951995   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#2016####
#visitations
m <- glm(rate.per.flower~net.treatment*insect.RTU + mean.temp + offset(mean.time), family = "quasipoisson", data = expt2.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df  Deviance Resid. Df Resid. Dev Pr(>Chi)  
## NULL                                         64    0.36189           
## net.treatment             1 0.0225113        63    0.33938  0.04503 *
## insect.RTU                1 0.0216441        62    0.31773  0.04937 *
## mean.temp                 1 0.0014249        61    0.31631  0.61408  
## net.treatment:insect.RTU  1 0.0285041        60    0.28780  0.02411 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#all contrasts
lsmeans(m, pairwise~net.treatment*insect.RTU, adjust="tukey")

## $lsmeans
##  net.treatment                  insect.RTU    lsmean        SE df
##  Larrea flowers with annuals    bees       -4.720968 0.1799804 NA
##  Larrea flowers without annuals bees       -4.077820 0.1652632 NA
##  Larrea flowers with annuals    flies      -4.674246 0.2899155 NA
##  Larrea flowers without annuals flies      -5.120969 0.3813389 NA
##  asymp.LCL asymp.UCL
##  -5.073723 -4.368213
##  -4.401730 -3.753910
##  -5.242470 -4.106022
##  -5.868380 -4.373559
## 
## Results are given on the log (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                                                                  
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,bees    
##  Larrea flowers with annuals,bees - Larrea flowers with annuals,flies      
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,flies   
##  Larrea flowers without annuals,bees - Larrea flowers with annuals,flies   
##  Larrea flowers without annuals,bees - Larrea flowers without annuals,flies
##  Larrea flowers with annuals,flies - Larrea flowers without annuals,flies  
##     estimate        SE df z.ratio p.value
##  -0.64314775 0.2856244 NA  -2.252  0.1096
##  -0.04672132 0.3055275 NA  -0.153  0.9987
##   0.40000139 0.4382498 NA   0.913  0.7981
##   0.59642643 0.3713445 NA   1.606  0.3751
##   1.04314914 0.3954240 NA   2.638  0.0415
##   0.44672272 0.4967443 NA   0.899  0.8052
## 
## Results are given on the log (not the response) scale. 
## P value adjustment: tukey method for comparing a family of 4 estimates

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -12.96, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##      alpha 
## -0.9949591

#require(multcomp)
#ph <- glht(m) #tests different from no effect
#summary(ph)

#rep at all as a nested, random effect
m <- glm(rate.per.flower~net.treatment*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "quasipoisson", data = expt2.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: quasipoisson, link: log
## 
## Response: rate.per.flower
## 
## Terms added sequentially (first to last)
## 
## 
##                              Df Deviance Resid. Df Resid. Dev Pr(>Chi)  
## NULL                                            64    0.36189           
## net.treatment                 1 0.022511        63    0.33938  0.03914 *
## mean.temp                     1 0.000228        62    0.33915  0.83542  
## insect.RTU:rep                7 0.070855        55    0.26829  0.06311 .
## net.treatment:insect.RTU:rep  6 0.018630        49    0.24966  0.74113  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#Overdisperion
m <- update(m, family = "poisson")
require(AER)
dispersiontest(m, trafo = 1) #not overdispersed

## 
##  Overdispersion test
## 
## data:  m
## z = -15.463, p-value = 1
## alternative hypothesis: true alpha is greater than 0
## sample estimates:
##      alpha 
## -0.9975496

#visit duration####
m <- glm(net.visitation~net.treatment*insect.RTU + mean.temp + offset(mean.time), family = "gaussian", data = expt2.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                          Df  Deviance Resid. Df Resid. Dev Pr(>Chi)   
## NULL                                         64   0.040293            
## net.treatment             1 0.0038412        63   0.036452 0.009179 **
## insect.RTU                1 0.0003985        62   0.036053 0.401397   
## mean.temp                 1 0.0002629        61   0.035790 0.495497   
## net.treatment:insect.RTU  1 0.0018358        60   0.033955 0.071688 . 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lsmeans(m, pairwise~net.treatment*insect.RTU, adjust="tukey")

## $lsmeans
##  net.treatment                  insect.RTU      lsmean          SE df
##  Larrea flowers with annuals    bees       0.004228414 0.005642967 NA
##  Larrea flowers without annuals bees       0.030781612 0.006292023 NA
##  Larrea flowers with annuals    flies      0.009411633 0.009283722 NA
##  Larrea flowers without annuals flies      0.010660084 0.009349877 NA
##     asymp.LCL  asymp.UCL
##  -0.006831598 0.01528843
##   0.018449474 0.04311375
##  -0.008784128 0.02760739
##  -0.007665339 0.02898551
## 
## Results are given on the identity (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                                                                  
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,bees    
##  Larrea flowers with annuals,bees - Larrea flowers with annuals,flies      
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,flies   
##  Larrea flowers without annuals,bees - Larrea flowers with annuals,flies   
##  Larrea flowers without annuals,bees - Larrea flowers without annuals,flies
##  Larrea flowers with annuals,flies - Larrea flowers without annuals,flies  
##     estimate          SE df z.ratio p.value
##  -0.02655320 0.009840541 NA  -2.698  0.0351
##  -0.00518322 0.009647604 NA  -0.537  0.9500
##  -0.00643167 0.011642336 NA  -0.552  0.9459
##   0.02136998 0.012539160 NA   1.704  0.3214
##   0.02012153 0.010319343 NA   1.950  0.2073
##  -0.00124845 0.013919974 NA  -0.090  0.9997
## 
## P value adjustment: tukey method for comparing a family of 4 estimates

#rep at all as a nested, random effect
m <- glm(net.visitation~net.treatment*insect.RTU%in%rep + mean.temp + offset(mean.time), family = "gaussian", data = expt2.2016)
anova(m, test = "Chisq")

## Analysis of Deviance Table
## 
## Model: gaussian, link: identity
## 
## Response: net.visitation
## 
## Terms added sequentially (first to last)
## 
## 
##                              Df  Deviance Resid. Df Resid. Dev Pr(>Chi)  
## NULL                                             64   0.040293           
## net.treatment                 1 0.0038412        63   0.036452  0.01322 *
## mean.temp                     1 0.0001690        62   0.036283  0.60327  
## insect.RTU:rep                7 0.0032357        55   0.033047  0.63902  
## net.treatment:insect.RTU:rep  6 0.0023906        49   0.030657  0.70087  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

lsmeans(m, pairwise~net.treatment*insect.RTU, adjust="tukey")

## $lsmeans
##  net.treatment                  insect.RTU      lsmean          SE df
##  Larrea flowers with annuals    bees       0.004397152 0.006032746 NA
##  Larrea flowers without annuals bees       0.031308301 0.006657298 NA
##  Larrea flowers with annuals    flies      0.003525286 0.011392158 NA
##  Larrea flowers without annuals flies               NA          NA NA
##     asymp.LCL  asymp.UCL
##  -0.007426813 0.01622112
##   0.018260236 0.04435637
##  -0.018802934 0.02585351
##            NA         NA
## 
## Results are averaged over the levels of: rep 
## Results are given on the identity (not the response) scale. 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                                                                  
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,bees    
##  Larrea flowers with annuals,bees - Larrea flowers with annuals,flies      
##  Larrea flowers with annuals,bees - Larrea flowers without annuals,flies   
##  Larrea flowers without annuals,bees - Larrea flowers with annuals,flies   
##  Larrea flowers without annuals,bees - Larrea flowers without annuals,flies
##  Larrea flowers with annuals,flies - Larrea flowers without annuals,flies  
##      estimate         SE df z.ratio p.value
##  -0.026911149 0.01049204 NA  -2.565  0.0505
##   0.000871866 0.01156498 NA   0.075  0.9998
##            NA         NA NA      NA      NA
##   0.027783015 0.01464250 NA   1.897  0.2291
##            NA         NA NA      NA      NA
##            NA         NA NA      NA      NA
## 
## Results are averaged over the levels of: rep 
## P value adjustment: tukey method for comparing a family of 4 estimates

#temp####
#bees only
fit <- lm(rate.per.flower~mean.temp, data = subset(expt1.2015, insect.RTU == "bees"))
summary(fit) 

## 
## Call:
## lm(formula = rate.per.flower ~ mean.temp, data = subset(expt1.2015, 
##     insect.RTU == "bees"))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.6799 -0.3619 -0.1014  0.1858  3.6439 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)   
## (Intercept) -1.09279    0.52235  -2.092  0.04067 * 
## mean.temp    0.05736    0.01942   2.953  0.00448 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6125 on 60 degrees of freedom
## Multiple R-squared:  0.1269, Adjusted R-squared:  0.1124 
## F-statistic: 8.722 on 1 and 60 DF,  p-value: 0.004484

fit <- lm(rate.per.flower~mean.temp, data = subset(expt1.2016, insect.RTU == "bees"))
summary(fit)

## 
## Call:
## lm(formula = rate.per.flower ~ mean.temp, data = subset(expt1.2016, 
##     insect.RTU == "bees"))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.43146 -0.77517 -0.04601  0.42583  3.06318 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)  
## (Intercept) -0.08360    0.69210  -0.121   0.9041  
## mean.temp    0.05632    0.02704   2.083   0.0404 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.006 on 82 degrees of freedom
## Multiple R-squared:  0.05024,    Adjusted R-squared:  0.03866 
## F-statistic: 4.338 on 1 and 82 DF,  p-value: 0.0404

#floral density####
#bees only
#expt1
fit <- lm(count~net.floral.density, data = subset(expt1.2015, insect.RTU == "bees"))
summary(fit)

## 
## Call:
## lm(formula = count ~ net.floral.density, data = subset(expt1.2015, 
##     insect.RTU == "bees"))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -13.667  -4.055  -1.611   2.650  27.246 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        4.743740   1.448461   3.275  0.00176 ** 
## net.floral.density 0.013923   0.002838   4.906 7.44e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 6.76 on 60 degrees of freedom
## Multiple R-squared:  0.2863, Adjusted R-squared:  0.2744 
## F-statistic: 24.07 on 1 and 60 DF,  p-value: 7.44e-06

a1 <- AIC(fit)

fit <- lm(count~net.floral.density, data = subset(expt1.2016, insect.RTU == "bees"))
summary(fit) 

## 
## Call:
## lm(formula = count ~ net.floral.density, data = subset(expt1.2016, 
##     insect.RTU == "bees"))
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -24.3698  -4.0359  -0.7851   3.9389  20.5886 
## 
## Coefficients:
##                    Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        5.356093   1.202305   4.455 2.63e-05 ***
## net.floral.density 0.036332   0.001472  24.677  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 7.598 on 82 degrees of freedom
## Multiple R-squared:  0.8813, Adjusted R-squared:  0.8799 
## F-statistic:   609 on 1 and 82 DF,  p-value: < 2.2e-16

a2 <- AIC(fit)

#expt2
fit <- lm(count~net.floral.density, data = subset(expt2.2015, insect.RTU == "bees"))
summary(fit)

## 
## Call:
## lm(formula = count ~ net.floral.density, data = subset(expt2.2015, 
##     insect.RTU == "bees"))
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -1.1216 -1.0136 -0.3366  0.0468  5.6429 
## 
## Coefficients:
##                     Estimate Std. Error t value Pr(>|t|)    
## (Intercept)        1.2295394  0.3534629   3.479 0.000963 ***
## net.floral.density 0.0044603  0.0006414   6.954  3.5e-09 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.611 on 58 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.4547, Adjusted R-squared:  0.4453 
## F-statistic: 48.36 on 1 and 58 DF,  p-value: 3.503e-09

a3 <- AIC(fit)

fit <- lm(count~net.floral.density, data = subset(expt2.2016, insect.RTU == "bees"))
summary(fit) 

## 
## Call:
## lm(formula = count ~ net.floral.density, data = subset(expt2.2016, 
##     insect.RTU == "bees"))
## 
## Residuals:
##        Min         1Q     Median         3Q        Max 
## -1.117e-14  2.607e-16  2.607e-16  2.969e-16  7.192e-16 
## 
## Coefficients:
##                      Estimate Std. Error    t value Pr(>|t|)    
## (Intercept)        -2.512e-15  3.850e-16 -6.525e+00 3.94e-08 ***
## net.floral.density  5.000e-03  6.034e-19  8.287e+15  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.692e-15 on 48 degrees of freedom
## Multiple R-squared:      1,  Adjusted R-squared:      1 
## F-statistic: 6.867e+31 on 1 and 48 DF,  p-value: < 2.2e-16

a4 <- AIC(fit)

#contrasts between annual and shrub flowers
#Quick interactive-term vs additive-term model comparison
x <- c(a1,a2)
x

## [1] 416.8747 583.0381

y <- c(a3, a4)
y

## [1]   231.4943 -3255.4590

z <- abs(x-y)
z

## [1]  185.3803 3838.4971

Interpretation

1. There is support for the single-magnet hypothesis through bees.
2. Both shrub species can generally facilitate annuals through pollination through bees.
3. Taxa specificity in plant-pollinator associations is critical for both foundation and annual plant species
4. Larrea flowers are less attractive to bees relative to annuals within this ecosystem.
   
5. Bees repond to temp by increasing visitation rates to annuals.
6. Bees respond to positively to increasing floral density of both annual and Larrea flowers.
7. There is evidence that shrubs receive no reciprocal benefit or facilitation from annuals, and furthermore, that there is a cost to facilitation by these foundation shrubs.
8. Bees are the most effective pollinators, i.e. they also spend the most time with annual flowers and likely are the most effective service pollinators.